% =========================================================================
% function priorstrv_loop
% Having obtained npoints draws, run csminwel for each of them
% Creates cell with all information
% Last Modified Januray 11  2008
% If MINTYPE == 5 use GENSYSPOST_FMINCON.M
% TOTCOF         dimension all coefficients
% NUMPAR/NPAR    dimension estimates ONLY
% =========================================================================

% Having the matrices
% PARSTRVALS
% LOGPRIOR_MAT
% LOGPOSTERIOR_MAT
% LOGLIKEL_MAT

% 1. Possibly start a DIARYFILE
if flag_diary==1
    diaryfname=['log ',strdate(1)];
    % Create a DIary File of the Optimization
    delete(diaryfname);
    diary(diaryfname);
    disp('__________________________________________________');
    datestr(now);
end
logprior_mat(:,2)= -1e15*ones(nstrvals,1);
loglikel_mat(:,2)= -1e15*ones(nstrvals,1);

% PARMODEMAT     Matrix of output vectors
% OPTIMAT        Matrix with info on optimization
%                First columns number of iterations
%                Second column return flags
%                Third time in minutes
parmodemat = zeros(numpar, nstrvals);
optimat = zeros(3, nstrvals);

% 2.
% Begin estimation Loop
% Will work with TRY and CATCH later on but not yet for debugging


% 3. Start a cell that can save output as we go along
% TAB_LOOP
%
tab_loop=emptycell(numpar+7,nstrvals+1);
tab_loop(1:numpar,1)=parnames;
tab_loop(numpar+2:end,1)={'logPost';'logLikel';'Number';'Iteration count';'exitflag';'minutes'};
cd(outpath);save tab_loop;cd(cucd);


% 4. Assign default values to ESTIMOPT if they do not exist
[e,estimopt]=ch_field(estimopt,'TolFun',1e-5);
[e,estimopt]=ch_field(estimopt,'MaxFuneval',200000);
[e,estimopt]=ch_field(estimopt,'MaxIter',20000);
[e,estimopt]=ch_field(estimopt,'Display','off');

%
% MaxFunIter
% Display
partrmat_est = partrmat(parposest, :);

% 5.
% PARSTRVALS(:,ii)
% ii=1;
% errcount=0;
for ii=1:nstrvals;
    disp('===================================================');
    disp('===================================================');
    disp(['Starting Candidate Value number : ',num2str(ii) ] );
    pause(0.5);
    tic;
    parstrval = parstrvals(:, ii);
    switch case_estima
        case 3      % CSMINWEL
            flag_transform = 1;
%             cd('O:\PROJ_LIB\LJM\exog\orig\lindu\april 2010 sim2'); 
%             load truepar fullp; 
%             x0mode=fullp(parposest)-0.01; 
%             
            x0mode = parstrval(parposest);
            
            % Recall former variable INDMAT must be trimmed
            % to ONLY estimated parameters
            x0min = mod2min(x0mode(:), partrmat_est); % MOD2MIN
            crit = estimopt.TolFun;
            nit = estimopt.MaxIter;
            nest=4*round(estimopt.MaxIter/4);
            dispaj('Running SA for ',nest,' iterations');
            
            [fh,xestmin,itct,retcode,fvec]=simulated_simp(@lmj_postoptim,x0min(:),5,1e5,nest,parstrval,parposest,partrmat_est,funcmod,Y,trainvec,prpar,prss,solveopt,addsol,ssposest,1);
            
            try
                figure;
                subplot(2,2,1);
                plot( fvec(1:round(nest/4)) );
                subplot(2,2,2);
                plot( fvec(round(nest/4):round(nest/2)));
                subplot(2,2,3);
                plot( fvec(round(nest/2):round(3*nest/4)));
                subplot(2,2,4);
                plot( fvec(round(3*nest/4):end) );
                
                cd(outpath);
                eval(['save fvec',num2str(ii),' fvec']);
                cd(cucd);
            catch
                
            end
            
 %eval( ['saveas(gcf,['Trace ',num2str(ii)],'fig'); 
 
            
            
            %                                         csminwel(', x0min, 0.5*eye(length(x0min)),[], crit, nit, ...
            %                             parstrval, parposest, partrmat_est, funcmod, Y, trainvec, prpar, prss, solveopt, ...
            %                             addsol, ssposest, flag_transform);
            
            
            
            %             [csmpost,xestmin,junk,H,itct,junk,retcode]=...
            %                 csminwel('lmj_postmin', x0min, 0.5*eye(length(x0min)),[], crit, nit, ...
            %                 parstrval, parposest, partrmat_est, funcmod, Y, trainvec, prpar, prss, solveopt, ...
            %                 addsol, ssposest, flag_transform);
            xestmode = min2mod(xestmin, partrmat_est);
            estimiter = itct;
            exitflag = retcode;
            
        case 2 % FMINCON
            
            flag_transform = 0;
            x0mode = parstrval(parposest);
            [xestmode, junk, exitflag, fminconout] = fmincon(@(x)lmj_postoptim(x,parstrval,parposest, ...
                partrmat_est, funcmod, Y, trainvec, prpar, prss, solveopt, addsol, ssposest, flag_transform), ...
                x0mode,[],[],[],[],prpar_lb(parposest),prpar_ub(parposest),[], estimopt);
            estimiter = fminconout.iterations;
            %         case 3
            %             x0mode = parstrval(parposest);            
            
        case 1
            
            flag_transform = 1;
            
            x0mode = parstrval(parposest);
            
            % Recall former variable INDMAT must be trimmed
            % to ONLY estimated parameters
            partrmat_est = partrmat(parposest, :);
            x0min = mod2min(x0mode(:), partrmat_est); % MOD2MIN
            crit = estimopt.TolFun;
            nit = estimopt.MaxIter;
            
            [csmpost,xestmin,junk,H,itct,junk,retcode]=...
                csminwel('lmj_postoptim', x0min, 0.25*eye(length(x0min)),[], crit, nit, ...
                parstrval, parposest, partrmat_est, funcmod, Y, trainvec, prpar, prss, solveopt, ...
                addsol, ssposest, flag_transform);
            
            xestmode = min2mod(xestmin, partrmat_est);
            estimiter = itct;
            exitflag = retcode;
            
    end
    parvecmode = zeros(numpar, 1);
    parvecmode(parposcal) = parcalval;
    parvecmode(parposest) = xestmode;
    
    parestmode = xestmode;
    
    
    [lpostd, likel] = lmj_postmax(parestmode, parvecmode, parposest, funcmod, Y, trainvec, prpar, prss, solveopt, addsol, ssposest);
    
    logpostd_mat(ii, 2) = lpostd;
    loglikel_mat(ii, 2) = likel;
    logprior_mat(ii, 2) = lpostd - likel;
    
    % PARMODEMAT     Matrix of output vectors
    parmodemat(:, ii) = parvecmode;
    
    % OPTIMAT        Matrix with info on optimization
    %                First columns number of iterations
    %                Second column return flags
    %
    
    % Print ESTIMATED Value and LogPost and LogLikelihood
    %
    %diary fmin
    
    time=toc/60;
    
    % Write it to optimat;
    optimat(1, ii) = estimiter;
    optimat(2, ii) = exitflag;
    optimat(3, ii) = time;
    
    % Write it to TAB_LOOP
    tab_loop(1:numpar,ii+1)=num2cprec(parmodemat(:,ii));
    tab_loop{numpar+2,ii+1}=num2str(logpostd_mat(ii,2));
    tab_loop{numpar+3,ii+1}=num2str(loglikel_mat(ii,2));
    tab_loop{numpar+4,ii+1}=num2str(ii);
    tab_loop{numpar+5,ii+1}= num2str(estimiter); % Iteration count
    tab_loop{numpar+6, ii+1}= num2str(exitflag); % Exit flag
    tab_loop{numpar+7, ii+1}= num2str(time);
    
    % ===============================
    % Save output of minimization
    % ===============================
    try
        cd(outpath);
        % SAVE % LOGPRIOR_MAT
        % LOGPOSTERIOR_MAT
        % LOGLIKEL_MAT
        save mat_loop logprior_mat logpostd_mat logprior_mat ...
            parmodemat optimat
        
        xlswrite('tab_loop', tab_loop);
        cd(cucd);
    catch
        disp('Could not save file')
    end
    
    % XLSWRITE
    % =====================================================================
    %    catch
    %         errcount=errcount+1;
    %         cd(outfold);
    %         diary errors;
    %         disp('Error in priorstrv_loop');
    %         dispaj('Iteration Number ',ii);
    %         lasterr
    %         diary off;
    %         load strv_all;
    %         postmodemat(ii)=-1e15;
    % =====================================================================
end
disp(' ');
disp('____________________________________');
disp('End all draws');
dispaj('Time in minutes=',toc/60);

diary off;
copyfile(diaryfname,outpath);
cd(cucd);